package pluginUtils;


public class Complexe {
	private double real,im,module,argument;
	
	public static final Complexe ZERO=createComplexeFromRI(0,0);
	public static final Complexe INFINITY=createComplexeFromRI(Double.POSITIVE_INFINITY,0); 
	public static final Complexe UN=createComplexeFromRI(1,0);
	public static final Complexe I=createComplexeFromRI(0,1);
	
	public static Complexe createComplexeFromRI(double r,double i){
		double mo=Math.sqrt(r*r+i*i); 
		double ar=Math.atan2(i,r); 	
		return new Complexe(r,i,mo,ar); 
	}
	
	public static Complexe createComplexeFromPolar(double m,double a){
		double r=m*Math.cos(a); 
		double i=m*Math.sin(a); 
		return new Complexe(r,i,m,a);
	}
	
	
	
	private Complexe(double real, double im, double module, double argument) {
		super();
		this.real = real;
		this.im = im;
		this.module = module;
		this.argument = argument;
	}
	
	public Complexe(Complexe source){
		this.real=source.real; 
		this.im=source.im; 
		this. module=source.module; 
		this.argument=source.argument; 
	}
	
	public boolean isInfinite(){
		return (Math.abs(this.module)==Double.POSITIVE_INFINITY); 	
	}
	
	public static Complexe sum(Complexe z1,Complexe z2){
		double r,i; 
		r=z1.real+z2.real; 
		i=z1.im+z2.im; 
		return createComplexeFromRI(r,i); 
	}
	
	public static Complexe sub(Complexe z1,Complexe z2){
		double r,i; 
		r=z1.real-z2.real; 
		i=z1.im-z2.im; 
		return createComplexeFromRI(r,i); 
	}
	
	public static Complexe mul(Complexe z1,Complexe z2){
		double r=z1.real*z2.real-z1.im*z2.im; 
		double i=z1.im*z2.real+z2.im*z1.real; 
		return createComplexeFromRI(r, i); 
	}
	
	public static Complexe opp(Complexe z1){
		return createComplexeFromRI(-z1.getReal(),-z1.getIm()); 
	}
	
	public static Complexe div(Complexe z1,Complexe z2){
		double r=z1.real*z2.real+z1.im*z2.im; 
		double i=z1.im*z2.real-z2.im*z1.real;
		double ratio=z2.squareNorm();
		return createComplexeFromRI(r/ratio, i/ratio); 
	}
	
	public static Complexe mulbyDouble(Complexe z,double x){
		return createComplexeFromRI(z.real*x,z.im*x);
	}
	
	public double squareNorm(){
			return this.real*this.real+this.im*this.im;
	}
	
	public double norme(){return Math.sqrt(squareNorm()); }
	
	public static Complexe log(Complexe z){
		return createComplexeFromRI(Math.log(z.module),z.argument); 
	}
	
	public static Complexe conjugate(Complexe z){
		return createComplexeFromRI(z.real, -z.im);
	}
	
	public Complexe conjugate(){
		return createComplexeFromRI(this.real, -this.im);
	}
	
	public static Complexe inverse(Complexe z){
		return Complexe.mulbyDouble(Complexe.conjugate(z), 1/z.squareNorm());
		
	}
	
	public static double distance(Complexe z1,Complexe z2){
		double d=(z1.real-z2.real)*(z1.real-z2.real)+(z1.im-z2.im)*(z1.im-z2.im); 
		return Math.sqrt(d); 
	}
	
	public static Complexe power(Complexe z1,Complexe z2){
		// TODO : z1^z2 ne marche pas quand z1=zero
		if((z1.equals(ZERO))&&!z2.equals(ZERO)) return ZERO; 
		double re,im; 
		double inner=z2.im*Math.log(z1.module)+z2.real*z1.argument;
		double coef=Math.pow(z1.module,z2.real)*Math.exp(-z2.im*z1.argument); 
		re=Math.cos(inner)*coef; 
		im=Math.sin(inner)*coef; 
		return createComplexeFromRI(re,im); 
	}
	
	public static Complexe sqrt(Complexe z){
		return createComplexeFromPolar(Math.sqrt(z.module),z.argument/2);
	}

	public double getReal() {
		return real;
	}

	public void setReal(double real) {
		this.real = real;
		this.module=Math.sqrt(this.real*this.real+this.im*this.im); 
		this.argument=Math.atan2(this.im,this.real); 	
	}

	public double getIm() {
		return im;
		
	}

	public void setIm(double im) {
		this.im = im;
		this.module=Math.sqrt(this.real*this.real+this.im*this.im); 
		this.argument=Math.atan2(this.im,this.real); 	
	}

	public double getModule() {
		return module;
	}

	public void setModule(double module) {
		this.module = module;
		this.real=this.module*Math.cos(this.argument); 
		this.im=this.module*Math.sin(this.argument); 
	}

	public double getArgument() {
		return argument;
	}

	public void setArgument(double argument) {
		this.argument = argument;
		this.real=this.module*Math.cos(this.argument); 
		this.im=this.module*Math.sin(this.argument); 
	}
	
	public String toString(){
			return "("+real+","+im+")-("+module+","+argument/Math.PI*180+")"; 
	}

	public int hashCode(){
		return new Double(this.real).hashCode()+(new Double(this.im).hashCode()); 
	}

	@Override
	public boolean equals(Object obj) {
		if (this == obj)
			return true;
		if (obj == null)
			return false;
		if (getClass() != obj.getClass())
			return false;
		final Complexe other = (Complexe) obj;
		if (Double.doubleToLongBits(im) != Double.doubleToLongBits(other.im))
			return false;
		if (Double.doubleToLongBits(real) != Double
				.doubleToLongBits(other.real))
			return false;
		return true;
	}
	
	/*
	from cmath import *
	 
	# Coefficients used by the GNU Scientific Library
	g = 7
	p = [0.99999999999980993, 676.5203681218851, -1259.1392167224028,
	     771.32342877765313, -176.61502916214059, 12.507343278686905,
	     -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7]
	 
	def gamma(z):
	    z = complex(z)
	    # Reflection formula
	    if z.real < 0.5:
	        return pi / (sin(pi*z)*gamma(1-z))
	    else:
	        z -= 1
	        x = p[0]
	        for i in range(1, g+2):
	            x += p[i]/(z+i)
	        t = z + g + 0.5
	        return sqrt(2*pi) * t**(z+0.5) * exp(-t) * x
*/
	
	public static Complexe sin(Complexe z){
		Complexe E=createComplexeFromRI(Math.E, 0); 
		Complexe one=power(E,mul(z,I));
		Complexe two=power(E,mul(z,opp(I)));
		return div(sub(one,two),mulbyDouble(I, 2)); 
		}
	
	public static Complexe cos(Complexe z){
		Complexe E=createComplexeFromRI(Math.E, 0); 
		Complexe one=power(E,mul(z,I));
		Complexe two=power(E,mul(z,opp(I)));
		return div(sum(one,two),createComplexeFromRI(2, 0)); 
		
	}
	
	public static Complexe gamma(Complexe z){
		int g=7; 
		double p[]={0.99999999999980993, 676.5203681218851, -1259.1392167224028,
		   	     771.32342877765313, -176.61502916214059, 12.507343278686905,
			     -0.13857109526572012, 9.9843695780195716e-6, 1.5056327351493116e-7};
		if(z.getReal()<0.5) {
				Complexe unmoinsz=sub(UN,z);
				return div(createComplexeFromRI(Math.PI,0),mul(sin(mulbyDouble(sub(UN,z),Math.PI)),gamma(sub(UN,z))));
		}
		
		z=sub(z,UN); 
		Complexe x=createComplexeFromRI(p[0],0);  
		for(int i=1;i<g+2;i++){
			Complexe coef=div(createComplexeFromRI(p[i],0),sum(z,createComplexeFromRI(i,0)));
			x=sum(x,coef); 
		}
			Complexe t=sum(z,createComplexeFromRI(g+0.5, 0));
			Complexe resu=createComplexeFromRI(Math.sqrt(2*Math.PI), 0); 
			resu=mul(resu,power(t,sum(z,createComplexeFromRI(0.5, 0)))); 
			resu=mul(resu,x);
			resu=mul(resu,power(createComplexeFromRI(Math.E, 0),opp(t)));
			return resu; 
	}
	
		private static double fact(int n){
			double sum=1.0; 
			for(int i=1;i<=n;i++)sum*=i; 
			return sum; 
		}

		public static Complexe airyA(Complexe z){
			int n=30; 
			double initCoef=1.0/(Math.PI*Math.pow(3,2.0/3.0));
			Complexe cumul=ZERO; 
			for(int i=0;i<n;i++){
				// calcul coef i
				double innerCoef=Math.sin(2.0*(i+1.0)*Math.PI/3.0); //ok
				Complexe terme=mulbyDouble(z,Math.pow(3.0,1.0/3.0));
				
				terme=power(terme,createComplexeFromRI(i+0.0,0));//ok
				
				terme=mul(terme,gamma(createComplexeFromRI(1.0/3.0*(i+1.0),0))); 
				
				terme=mulbyDouble(terme,1.0/fact(i)*innerCoef);
				
				cumul=sum(cumul,terme); 
				
			}
			return mulbyDouble(cumul,initCoef); 
		}
		
		public static Complexe airyB(Complexe z){
			Complexe w=createComplexeFromPolar(1, 5*Math.PI/7.0); 
			Complexe w2=mul(w,w); 
			Complexe q=createComplexeFromRI(7, 0); 
			
			Complexe q2=mul(q,q); 
			Complexe terme1=mul(I,q2); 
			terme1=mul(terme1,airyA(mul(z,w2)));
			Complexe terme2=mulbyDouble(mul(q,I),-1);
			terme2=mul(terme2,airyA(mul(z,w))); 
			return sum(terme1,terme2);
		}

		public static Complexe powerDouble(Complexe z, double exponent) {
			if(Math.abs(exponent)<1e-10) return Complexe.UN; 
			Complexe z1=Complexe.createComplexeFromRI(exponent,0); 
			return Complexe.power(z,z1);
			
		}
	
}
